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Due to remarkable advances in colloid synthesis techniques, systems of squares and cubes, once 
an academic abstraction for theorists and simulators, are nowadays an experimental reality. By 
means of a free minimization of the free-energy functional, we apply Fundamental Measure Theory 
to analyze the phase behavior of parallel hard squares and hard cubes. We compare our results with 
those obtained by the traditional approach based on the Gaussian parameterization, finding small 
deviations and good overall agreement between the two methods. For hard squares our predictions 
feature at intermediate packing fraction a smectic phase, which is however expected to be unstable 
due to thermal fluctuations. This implies that for hard squares the theory predicts either a vacancy- 
rich second-order transition or a vacancy-poor weakly first-order phase transition at higher density. 
In accordance with previous studies, a second-order transition with a high vacancy concentration is 
predicted for hard cubes. 



I. INTRODUCTION 

Hard spheres represent the simplest and most versa- 
tile model for the description of molecular and colloidal 
many-particles systems. This statement is particularly 
true since 1957, when Wood and Jacobson ^1^ and Alder 
and Wainwright |2] demonstrated that hard spheres un- 
dergo a fluid-to-crystal transition, and therefore that 
hard interactions alone can account for freezing. 

Systems of hard cubes, on the other hand, were consid- 
ered as mere toy models until only a few years ago. The 
reason for this is evident: no molecule or macromolec- 
ular aggregate found in nature is known to be reason- 
ably approximated by this shape. However, the interac- 
tion between parallel hard cubes is the second-simplest 
hard interaction one can imagine after that between hard 
spheres. Its simplicity made this model a perfect object 
of study for theory and simulation. 

Early studies on the equation of state of parallel hard 
squares [D — 2 dimensions) and cubes {D ~ 3) date 
back to the dawn of computer simulation in the 1950s [3l 
13]. Soon after, the question regarding the high-packing 
phase behavior of the models arose. For parallel hard 
squares, a transition from the fluid to a square-lattice 
crystal (with quasi- long-range order) was found [5l|6], but 
its character, whether continuous or discontinuous, has 
been a matter of debate ever since [TlIIO). Conversely, 
the stability of a "brick-wall" smectic phase with one- 
dimensional ordering in rows (or columns) was suggested 
to exist for parallel hard squares, but the stability of this 
peculiar state was soon ruled out [TT] . Similarly, parallel 
hard cubes manifest a fluid-to-crystal transition with a 
well-established second-order character [12] and no stable 
phase with lower translational symmetry than the crystal 
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is expected [13]. 

In the mid-1990s, while the interest of the liquid-state 
community was focusing on mixtures of hard spheres, 
hard cubes were rediscovered. By means of computer 
simulation Dijkstra et al. showed evidence of a demix- 
ing transition in a binary system of parallel hard cubes 
on a lattice, thus demonstrating that additive hard in- 
teractions can induce an entropy-driven fluid-fluid phase 
separation [TH |TS] . These results motivated Cuesta and 
Martinez- Raton to face the problem by means of density 
functional theory. Following the pioneering approach de- 
veloped by Rosenfeld for hard spheres [Hj, they devel- 
oped a Fundamental Measure Theory (FMT) [T^ formal- 
ism aimed at describing both the homogeneous and in- 
homogeneous phase behavior of mixtures of squares and 
cubes [THHn]. 

Since the early work on hard squares and cubes, the 
progress in colloidal particles synthesis has been enor- 
mous. In particular, colloidal suspensions of micron-sized 
cubes [20] and quasi-two-dimensional square platelets 
[211 have been recently produced and analyzed. These 
experimental advances led to a renewed interest in the 
model, and at present more complex aspects like the role 
of orientational degrees of freedom, the addition of dipo- 
lar interactions, the roundedness of the shape and the 
effect of vacancies in the freezing mechanism constitute 
objects of intense research [22H27] . Far from being a toy 
model or a mere academic exercise, squares and cubes 
have therefore gained a key role as model systems of non- 
spherical colloidal particles. 

Besides the development of new theories [351 HH] , the 
increasing attention towards the self-assembly of non- 
spherical particles requires a detailed analysis of the ca- 
pabilities of the existing ones. The aim of this paper 
is to reinvestigate the prediction of Fundamental Mea- 
sure Theory as formulated in Ref. (I8j for the phase 
behavior of parallel hard squares and cubes. The focus 
of our attention points to the freezing transition and the 
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structure of the high-density inhomogeneous phases. In 
particular, by exploiting present-day computer power we 
improve previous analyses on the subject by performing a 
free minimization of the density functional, and compare 
our results with those obtained by means of the widely 
applied Gaussian parameterization of the single-particle 
density. We observe good overall agreement between the 
two methods, the main drawbacks of the Gaussian pa- 
rameterization being (i) a systematic albeit small under- 
estimation of the equilibrium vacancy concentration in 
the crystal and (ii) the lack of anisotropy of the crystal 
density peaks at high enough density. Furthermore, to 
the best of our knowledge this work constitutes the first 
density-functional theory study of the hard-square sys- 
tem. We show that Fundamental Measure Theory sur- 
prisingly predicts a smectic phase absent in computer 
simulation and suggests that, in analogy with the hard- 
cube system, vacancies can play a crucial role in the freez- 
ing transition. 



where {na{r)} is a set of weighted densities, labeled by a, 
obtained as convolutions between the single-particle den- 
sity and a set of corresponding weight functions Wjj(r), 



n„(r)- J d^rV(r')w„(r-r'). (3) 

The functional dependence of $''^^({^0}) is determined 
by extrapolating from known limiting cases, such as 
the homogeneous bulk equation of state, the low-density 
second-virial behavior, and the dimensional crossover to 
highly confined conditions |17) . 

For hard parallel squares {D = 2) and cubes {D = 
3) with side a, the FMT functional was determined by 
Cuesta and Martinez- Raton in Ref. [18] . In what follows, 
we report the explicit expression of ^^^\{na}) for the 
single-component case. Following Ref. [TH] we introduce 
the auxiliary functions 



II. DENSITY FUNCTIONAL THEORY 

The density-functional theory route to the equilibrium 
properties of a many-body system consists of expressing 
the intrinsic Helmholtz free energy as a functional of 
the single-particle density p(r) [30] . When considering a 
system composed of a single species of particles having 
only translational (and no rotational) degrees of freedom 
in D dimensions, the free-energy functional reads 

= J d^rp(r){log[p(r)A3] - 1 } + /3 [p] , (1) 

where r is a I?-dimensional vector, /3 — {ksT)'^ is the in- 
verse temperature in units of the Boltzmann constant, A 
the thermal wavelength and the integrals are performed 
over the (D-dimensional) volume V occupied by the sys- 
tem. The first term in the right-hand side of Eq. ([T]) 
denotes the ideal-gas contribution, while the second de- 
scribes the excess contribution due to particle-particle 
interactions. 



T(2:)-e(^|-|.T|), C(a:) = i5(^|-N), (4) 

defined for a; e M. 

For parallel squares the weight functions are 

wo(r) = C(a^)C(y); (5a) 
wi(r) = (C(x)r(y),T(x)C(y)); (5b) 
W2(r) = r(.T)T(2/), (5c) 

where we note that wi has a vector character. The func- 
tional dependence of the excess free energy of parallel 
cubes is given by 



$(2) = -nolog(l-,i2) + ^^^. (6) 

1-712 

For parallel cubes the weight functions are 



A. Fundamental Measure Theory (FMT) 

The excess free-energy functional J^^^'^ in Eq. ([T]) is 
the non-trivial element of the theory: it contains the free- 
energy dependence on the inter-particle interactions and 
it can not be calculated exactly in general. 

Various methods to systematically estimate this func- 
tional dependence have been developed. For hard spheres 
the undoubtedly most successful approach is that of Fun- 
damental Measure Theory (FMT). According to FMT, and 
the excess free energy is written as 



wo(r) 


= C(a:)C(y)C(^); 


(7a) 


wi(r) 


- (T(x)C(y)C(^),C(x)r(y)C(z; 


l,C(a;)C(y)r(z)); 






(7b) 


W2(r) 


{C{x)T{y)T{z),T{x)C{y)T{z 


),T(a;)T(y)C(z)); 






(7c) 


W3(r) 


= T{x)T{y)T{z); 


(7d) 



:iD. 



(x) iv) (z) 



3 



B. Functional minimization 

Once an explicit expression for the functional depen- 
dence of J^'^^'^ on the single-particle density p(r) is estab- 
lished, the equilibrium Helmholtz free energy F{T, V, N) 
of particles at temperature T in a volume V (area A 
for D = 2) is obtained as the minimum of Eq. ([T]) with 
respect to p{r) under the constraint that 



d^rp(r) = N. 



(9) 



In the case of hard-core interactions between the par- 
ticles, the system is athermal and its thermodynamic 
state is completely identified by the packing fraction 
1] = Nvp/V, where Vp = is the particle's volume (for 
D — 2 dimensions rj — Nap/ A and Op = cr^). 

The numerically easiest way to solve the functional 
minimization problem consists of expressing the single- 
particle density in terms of a limited number of varia- 
tional parameters. After inserting this ansatz into the 
free energy, the latter is minimized with respect to the 
variational parameters to obtain an estimate of the free 
energy at equilibrium. This approach has been widely ap- 
plied in studying the freezing transition of hard spheres, 
where the single-particle density was parameterized as a 
sum of Gaussian functions centered on the lattice sites of 
the expected stable crystal phase (Gaussian parameter- 
ization or ansatz) In this paper we investigate 
the freezing transition of squares and cubes into square 
and simple-cubic crystal phases, for which the Gaussian 
ansatz can be expressed as 



directions. In practice, these spatial modulations must 
be inserted explicitly into Eq. ( |11[ ) by means of a Fourier 
series expansion. Therefore, the single particle density 
of a phase characterized by a d-dimensional spontaneous 
breaking of the translational symmetry is obtained by 
solving the following equation 




exp 



5/3 J"^ 



Spis') 



(12) 

where s G F is a d-dimensional vector, F — [—A/2, A/2]'' 
and A is the periodicity of the inhomogeneous solution 
(assumed to be the same along all the d directions). The 
minimization procedure consists of (i) solving Eq. ( 12 ) at 
fixed A and (ii) identifying the value of A that minimizes 
the free energy Eq. For hard squares (D = 2) we 

will see that Eq. (12) describes smectic {Sm, d — 1) 
and square crystal JX, d = 2) phases; for hard cubes 
[D = 3) Eq. (12) accounts for smectic {Sm, d — 1), 
columnar [Col, d = 2) and simple-cubic crystal [X , d ~ 
3) ordering. 



A numerical solution of Eq. ( 12 ) on a grid of points is 



expected to offer a better description of the single particle 
density, and therefore a lower minimum free energy, than 
by minimizing the free energy by means of the Gaussian 
ansatz Eq. (10). We develop a Picard algorithm to solve 



Eq. ( 12 ), where all the convolutions involved in the FMT 



formalism are handled by means of Fast Fourier Trans- 
forms [36]. Moreover, the minimization with respect to 
the lattice spacing A is performed using the Brent algo- 
rithm j37j . 



In Eq. 



A\-D 



E 



exp 



-7(r - An)2 



(10) 



( 10 ) the variational parameters are half of the 
inverse variance 7 and the lattice constant A. Note that 
the lattice constant A is related to the vacancy concen- 
tration of the crystal Xyac = {Ns\tcs - N)/Nsitcs through 
Xvac = 1 — vi^/'^)^ 1 where A'sitcs is the total number of 
sites. 

An alternative to the Gaussian parameterization con- 
sists of the numerical solution of the Euler-Lagrange 
equation associated with the minimization problem (free 
minimization) \34\ 1351 . By performing a functional dif- 
ferentiation of Eq. (n]) and imposing the constraint of 
Eq. ([9]), one finds that the equilibrium p{r) satisfies the 
following self-consistency equation 



p(r) = N exp 



(11) 

At high packing fraction 77 one expects the free energy 
to be minimized by inhomogeneous solutions character- 
ized by spatial modulations of p{r) along one or more 



r 5/37"°'^'= 1 




r Sf3T^'''^i 


I 6p{r) \ 


|y d^r' exp 


[ 5p{r') \ 



III. RESULTS 

A. Parallel Hard Squares (D=2) 

When considering the high-density phase behavior, 
monodisperse squares (as well as cubes in _D = 3 dimen- 
sions) possess a peculiar property. Unlike other regular 
polygons (e.g., regular triangles, pentagons, hexagons,..., 
and disks), squares do not have a well-defined "locked- 
in" configuration at close packing. In other words, be- 
sides the two-dimensional ordered square crystal {X), 
any other configuration with rows (or columns) shifted 
with respect to one another completely fills the plane. 
Therefore, also a smectic phase (Sm), characterized by 
positional ordering along one direction only, should in 
principle be considered as a candidate stable phase (see 
Fig. [ija)). The higher degeneracy of 5*7?^ configurations 
with respect to X configurations suggests a higher en- 
tropy of the former with respect to the latter. How- 
ever, in low-dimensional systems thermal fluctuations 
from equilibrium can play a relevant role in destroying 
long-range order, leading to so-called Landau-Peierls in- 
stabilities |39j . In particular, for short-range interactions 
proper crystals do not exist in Z? = 2 dimensions, since 
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FIG. 2. FMT results for (a) the vacancy concentration and 
(b) the root-mean-squared deviation from the average posi- 
tion (in units of the lattice constant A) of the square crys- 
tal phase of parallel hard squares. Solid lines correspond 
to values calculated by free minimization of the FMT func- 
tional, whereas dashed lines indicate those obtained through 
the Gaussian ansatz. 




the single-particle density of the X phase by assuming 
long-range order (cf. Eqs. (10 1 and (12|). A similar 
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FIG. 1. (a) Pictorial representation of (from left to right) 
the fluid, smectic (Sm) and square crystal (X) phases of par- 
allel hard squares, (b) Phase diagram of parallel hard squares 
according to FMT, to be compared with the simulation value 
of Ref. [To] for the fluid-to-crystal transition packing frac- 
tion (vertical arrow). The (Sm) interval highlights the states 
where the (Peierls-Landau unstable) smectic phase is pre- 
dicted to be the stable phase, (c) FMT prediction for the 
free energy per unity area F/A of fluid (dotted green), Sm 
(dashed blue) and X (solid red lines) phases in the proximity 
of the second-order fluid-to-smectic (black circle) and (d) the 
first-order smectic-to-crystal transition. 



positional ordering can in this case have at most quasi- 
long-range character [JD]. The situation is even more 
dramatic when considering smectic phases, where ther- 
mal fluctuations make the correlation between layers de- 
cay exponentially with the distance [H] . This means that 
in D = 2 dimensions we do not expect smectic ordering 
to be stable at all. Computer simulations of both parallel 
[8] and freely-rotating [38] hard squares, where only a di- 
rect fluid-to-crystal phase transition was observed with- 
out any smectic state, confirm this picture. 

We report in Fig. [TJb) the phase diagram of parallel 
hard squares, as obtained by freely minimizing the FMT 
functional with respect to the single-particle density p(r) 
including vacancies. Despite the above-mentioned con- 
siderations on the effect of fluctuations, we approximate 



approximation was recently applied for the description 
of the freezing transition in two-dimensional hard disks, 
showing remarkably good agreement with computer sim- 
ulations [3S]. We include in our calculations also the 
possibility of long-range Sm ordering, which is however 
expected to be Landau-Peierls unstable. The free-energy 
dependence of the fluid, smectic and crystal phases on the 
packing fraction is reported in Fig. [ijc) and (d) for two 
density intervals. Note that this representation allows 
for common tangent construction to identify coexisting 
states. Surprisingly, FMT predicts a second-order fluid- 
to-smectic transition at 77* = 0.538 (Fig. [ijc)) and a 
weakly first-order smectic-to-crystal transition with bulk 
coexisting densities rjsm = 0.726 and r]x = 0-730 (Fig. 
[ijd)). The picture does not change appreciably by mini- 
mizing the free energy within the Gaussian ansatz, giving 
the sole effect of slightly displacing the alleged Sm — X 
transition (rjsm — 0.750 and rjx — 0.756, not shown). As 
already pointed out, theoretical considerations and sim- 
ulation results rule out the possibility of stable smectic 
ordering in the thermodynamic limit. Therefore, we must 
conclude that the smectic phase is an artifact due to the 
mean-field character of the Fundamental Measure The- 
ory, which is unable to take fully into account the role of 
long-wavelength fluctuations from equilibrium. On the 
other hand, the question whether the fluid or the crys- 
tal is the stable phase in the range of allegedly smectic 
stability (striped region in Fig. [ijb)) is open. We ad- 
dress this point, as well as possible conditions of smectic 
stability, in the final discussion of Sec. |IV[ 

In order to further investigate the properties of the 
crystal, we report the dependence of the vacancy concen- 
tration (Fig. [2|^a)) and the root-mean-squared deviation 
from the average position in the unit cell, also known as 
Lindemann parameter (Fig. [2|^b)), on the packing frac- 
tion. At the second-order transition at 77 = 0.538 the 
vacancy concentration is x^ac — 15%, 
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FIG. 3. Single-particle density p{x, y) in the unit cell of the 
square crystal phase obtained through FMT free minimization 
at packing fraction (a) r) = 0.55, (c) r] = 0.65 and (e) r) = 
0.75. On the right column we report the section of the single- 
particle density along the [10] (red lines) and [11] (blue lines) 
crystallographic directions. The graphs are expressed in log- 
scale as a function of the squared distance from the center of 
the unit cell. The dashed lines represent the corresponding 
functional dependence obtained through the Gaussian ansatz. 



confined regime at higher density. In Fig. [sjja), (c) and 
(e) we report the functional dependence of the equilib- 
rium single-particle density inside the unit cell at packing 
fraction -q — 0.55, 0.65 and 0.75, respectively. In order to 
ease the analysis, we plot on the right of each figure (Fig. 
[Sj^b), (d) and (f)) a section of the corresponding p{x,y) 
along the crystallographic directions [10] and [11]. These 
graphs are represented in logarithmic scale as a function 
of to highlight Gaussian behavior (straight lines). As 
expected, the section along the [10] direction, connect- 
ing nearest-neighbor sites, is systematically bigger than 
that along the [11] direction for both the freely mini- 
mized (solid lines) and Gaussian-parameterized (dashed 
lines) profiles. At the lowest packing fraction the peak of 
the density distribution is smeared out on the unit cell. 
As a consequence, the single-particle density in the Gaus- 
sian parameterization shows relevant deviations from the 
Gaussian distribution due to the overlap of the peaks cen- 
tered on neighboring cells. In this way, the tails of the 
neighboring lattice sites allow to properly account for the 
anisotropy of the density peak, thus leading to a marked 
difference between the [10] and [11] profiles, similar to 
the case of free minimization (Fig. Sj^b)). This overlap 
is weaker at higher packing fraction, where the confine- 
ment is stronger; hence, the Gaussian ansatz fails to re- 
produce the anisotropy of the distribution in this regime 
(Fig. [sjjd) and (f)). However, these deviations occur on a 
density scale a few orders of magnitude smaller than the 
peak value, and therefore their relevance is quantitatively 
limited. 



B. Parallel Hard Cubes (D=3) 



monotonically with rj. This value is appreciably higher 
than that predicted for hard disks Xyac — 2% [3S], in 
analogy to the D = 3 dimensional case of hard-cubes 
and hard-spheres systems, where the vacancy concentra- 
tion at the fiuid-coexistence in the former is two orders of 
magnitude higher than in the latter [HI [571 1311 • un- 
like the case of hard cubes [22] , these results suggest that 
vacancies can play an important role in stabilizing (quasi- 
long-range) crystal order in systems of hard squares. 
Moreover, Fig. [2ja) highlights a systematic, albeit small, 
underestimation of the vacancy concentration at inter- 
mediate packing fraction when the Gaussian ansatz is 
applied. On the other hand, no appreciable difference 
with the root-mean-squared deviation calculated by free 
minimization is observed in Fig. [2jb). In both cases, 
at the second-order transition the Lindemann parameter 
assumes a value close to 0.4, remarkably higher than that 
of 0.15 expected at melting for three-dimensional systems 
[42l . and decays to zero towards close packing. 

We complete our analysis of the crystal phase by study- 
ing the evolution of the single-particle density from the 
freezing-transition region, where p{x,y) is still apprecia- 
bly non-zero at the edge of the Wigner-Seitz cell, to the 



Here we compare the predictions of FMT for parallel 
hard cubes, extensively studied in the past by means of 
the Gaussian parameterization |13l [HI [27l , with our re- 
sults based on the free minimization of the functional. 
Since the formulation by Cuesta and Martfnez-Raton 
[18) . FMT is known to correctly predict two significant 
properties of the freezing transition of the model [19] : (i) 
its second-order character and (ii) the role of vacancies 
in stabilizing the crystal. 

The second-order fluid-to-crystal transition, which is 
known to become first-order when the rotational degrees 
of freedom are taken into account [T^l, is predicted to 
occur at 77 = 0.314. As in the case of parallel squares, 
this value appreciably underestimates the simulation re- 
sult of 77 = 0.469 ]57] (see Fig. |4]^a)). Also in analogy 
with parallel hard squares, parallel cubes lack a "locked- 
in" configuration at close packing. By means of a bi- 
furcation analysis of the FMT functional and computer 
simulations, Groh and Mulder addressed the question 
about the stability of columnar order, and showed it to 
be metastable. In contrast to hard squares, smectic and 
columnar solutions are in the three-dimensional case al- 
ways metastable with respect to the crystal. This finding, 
which results directly from our free minimization scheme. 
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FIG. 4. (a) Phase diagram of parallel hard cubes according 
to FMT, to be compared with the simulation result of Ref. 
|27| for the fluid-to-crystal (X) transition packing fraction 
(vertical arrow), (b) FMT prediction for the free energy per 
unity volume PF/V of fluid (dotted green), crystal (solid red) 
and the metastable columnar [Col, dashed blue) and smec- 
tic [Sm, dot-dashed orange) phases in the proximity of the 
second-order freezing transition (black circle). 





0.3 0.4 0.5 0.6 0.7 0.8 
11 



0.3 0.4 0.5 0.6 0.7 0.8 0.9 
11 



FIG. 5. FMT results for (a) the vacancy concentration and 
(b) the root-mean-squared deviation from the average posi- 
tion (in units of the lattice constant A) of the simple-cubic 
crystal phase of parallel hard cubes. Solid lines correspond 
to values calculated by free minimization of the FMT func- 
tional, whereas dashed lines indicate those obtained through 
the Gaussian ansatz. 



is easily verified by comparing the free-energy curves of 
smectic, columnar and crystal phases as a function of the 
packing fraction in Fig. |4f^b). 

The remarkably high concentration of vacancies at the 
freezing transition, x^ac — 30%, is a known feature of the 
theory (cf. Fig. [5ja)) [131 [T^. Despite this value being 
three orders of magnitude higher than that measured for 
hard spheres [34], it was shown to be compatible with 
computer simulations of both parallel [x^ac = 13% [27] ) 
and freely-rotating {x^ac = 6.4% P5]) hard cubes, thus 
highlighting the essential role of vacancies in stabilizing 
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FIG. 6. FMT prediction for the single-particle density of 
parallel hard cubes. The graphs on the left show sections 
of p{x, y, z) along the crystallographic directions [100] (red 
lines), [110] (green lines) and [111] (blue lines) at packing 
fraction (a) 77 = 0.32, (c) 0.50 and (e) 0.70 calculated by 
means of free minimization of the functional (solid lines) and 
the Gaussian ansatz (dashed lines). The right graphs repre- 
sent the corresponding absolute difference between the free 
minimization solution and the Gaussian ansatz for the three 
packing fractions ((b)?7 = 0.32, (d) 0.50 and (f) 0.70) and the 
three crystallographic directions considered. 



the simple-cubic crystal. Within the free minimization of 
the FMT functional, the vacancy concentration at bulk 
coexistence does not change with respect to the Gaussian 
ansatz result. Nonetheless, an inspection of Fig. [sja), 
reporting x^ac as a function of the packing fraction 77, 
shows that the Gaussian ansatz tends to underestimate 
this property at intermediate packing fractions. There- 
fore, in this regime the free minimization improves the 
Gaussian ansatz data by furnishing results closer to those 
of computer simulation [27^ . However, if we focus on the 
root-mean-squared deviation from the average lattice site 
(Lindemann parameter), reported as a function of rj in 
Fig. [5|b), we do not observe any appreciable deviation 
from the known dependence calculated by means of the 
Gaussian ansatz. 

Finally, in Fig. [6] we represent sections of the equi- 
librium single-particle density p{x, y, z) at packing frac- 
tion 77 = 0.32 ((a)-(b)), 0.50 ((c)-(d)) and 0.70 ((e)-(f)) 
and we compare them with the corresponding Gaussian 
ansatz solution (dashed lines). The graphs on the left 
((a), (c) and (e)) show sections along the crystallographic 
directions [100] (red lines), [110] (green lines) and [111] 



7 



(blue lines); to ease the comparison, we report on the 
right of each graph ((b), (d) and (f)) the absolute dif- 
ference of these sections between the two minimization 
methods. In the three cases, the cubic symmetry of 
the freely-minimized solution is evident by the hierarchy 
in values of the single particle density along the three 
crystallographic directions. In analogy with the parallel 
square system, at low enough packing fraction there is 
good quantitative agreement between the two methods, 
since the overlap between neighboring peaks within the 
Gaussian ansatz allows to reproduce the anisotropy of 
the single-particle distribution. At higher packing frac- 
tion, deviations from the Gaussian-ansatz solution are 
evident, but limited to the low-density region of the unit 
ceU. 



IV. DISCUSSION AND CONCLUSIONS 

By means of Fundamental Measure Theory we investi- 
gate the phase behavior of single-component systems of 
parallel hard squares (in D ~ 2 dimensions) and cubes 
{D — 3). Our attention focuses on the predictions for 
the freezing transition and the properties of the crystal 
phase. In density-functional theory the typical approach 
for describing crystal phases is based on the parameter- 
ization of the single-particle density by a sum of Gaus- 
sian functions centered on the lattice sites. We compare 
these predictions with a more accurate free-minimization 
method, where the single-particle density is evaluated on 
a grid of points. 

Despite its simplicity, we conclude that for both 
squares and cubes the Gaussian parameterization works 
remarkably well. Apart from some inadequacy of the 
Gaussian ansatz in describing the anisotropy of the 
single-particle density of the crystal, the main deviations 
between the two minimization methods lie in the ex- 
pected vacancy concentrations of the square and simple- 
cubic crystals, which appears to be slightly underesti- 
mated by the Gaussian ansatz. On the other hand, as 
already noticed for cubes, FMT suffers from a serious in- 
ability to give quantitatively reliable values for the freez- 
ing packing fraction. However, improvement in this di- 
rection can be achieved only by a reformulation of the 
theory itself, as the numerical minimization is performed 
exactly. 

For the three-dimensional system of parallel hard cubes 
our results coincide with previous FMT analysis based 
on the Gaussian parameterization and indicate a second- 
order vacancy-rich fluid-to-crystal transition. For the 



parallel hard-square system, this work constitutes to the 
best of our knowledge the first analysis based on density- 
functional theory. In contrast with previous simulation 
studies. Fundamental Measure Theory predicts a stable 
smectic phase in between the low-density fluid and the 
high-density square crystal. However, by taking into ac- 
count the effect of long-wavelength thermal fluctuations, 
one can show the one-dimensional smectic ordering to 
be unstable. Therefore, we deduce that the mean-field 
character of the theory, which is unable to properly take 
into account the role of fiuctuations from equilibirum, is 
the element to be blamed for this artifact. When big 
enough simulation boxes are considered, computer sim- 
ulations with periodic boundary conditions comfirm the 
picture of an unstable smectic phase f5F|. However, it 
is interesting to notice that, when the simulation box is 
small enough, ordering of the squares in parallel layers 
was observed [38]. On the basis of these observations, 
it is tempting to conclude that, when long- wavelength 
fiuctuations can be neglected, the behavior of the system 
coincides with the predictions of FMT, showing a stable 
smectic phase. In other words, we expect parallel-squares 
systems to develop intermediate smectic states in finite 
size systems and under the effect of confining walls. 

For what concerns the original problem regarding the 
phase behavior of parallel hard squares in the thermo- 
dynamic limit, the conclusions we can draw are more 
limited. In fact, on the basis of our theoretical results 
we do not have enough elements to deduce which of the 
two phases, either the fiuid or the crystal, is the sta- 
ble one in the density range where the theory predicts a 
stable smectic. Therefore, we cannot conclude whether 
the theory predicts a low-density second-order freezing 
transition with a high vacancy concentration in the crys- 
tal, or a higher-density weakly first-order freezing with a 
lower vacancy concentration. Addressing this problem is 
challenging, as it should involve a proper incorporation 
of long-wavelength Landau-Peierls fluctuations into the 
density-functional theory. 
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